rsyrs.est=function(data,yi,k,N,param){
  data=as.name(data); N=N; k=k
  ns=nrow(yi); n=k*ns
  fpc=(1-n/N); param=as.name(param)
  ybari=apply(yi[,-1],1,mean,na.rm=T)
  muhat=sum(ybari)/ns
  s2ybar=sum((ybari-muhat)^2)/(ns-1)
  if(param=='mean'){
    est=muhat
    vhat=(s2ybar/ns)*fpc
    bound=round(2*sqrt(vhat),2)
    lower=est-bound; upper=est+bound
  } else if(param=='total'){
    est=N*muhat
    vhat=N^2*(s2ybar/ns)*fpc
    bound=round(2*sqrt(vhat),2)
    lower=est-bound; upper=est+bound
  } 
  cat("","\n","Results from Repeated 1-in-",k,"SyRS: Data =",data,"\n","N =",N,"n =",n,"\n",
      "FPC =",fpc,"\n","Estimate of",param," =",est,"\n","Vhat =",vhat,"\n","Bound =",
      bound,"\n","Lower Bound =",lower,"Upper Bound =",upper,"\n","")
  results=list(est=est,data=data,n=n,ns=ns,N=N,fpc=fpc,vhat=vhat,param=param,bound=bound,lower=lower,upper=upper,s2ybar=s2ybar)
}
# to use the function with its call:
# rsyrs.est(data,yi,k,N,param)

# data: name of dataset, in quotes
# yi: a numeric vector or matrix containing the values (each row is one sample); first column is an index of the record number
# k: 1-in-k value 
# N: population size
# param: c('mean','total')